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Recent experimental results have revealed a surprisingly rich structure of multiparticle azimuthal 
correlations in high energy proton-nucleus collisions. Final state collective effects can be responsible 
for many of the observed effects, but it has recently been argued that a part of these correlations are 
present already in the wavefunctions of the colliding particles. We evaluate the momentum space 
2-particle cumulant azimuthal anisotropy coefficients v„{2},n = 2,3,4 from fundamental represen¬ 
tation Wilson line distributions describing the high energy nucleus. These would correspond to the 
flow coefficients in very forward proton-nucleus scattering. We find significant differences beteen 
Wilson lines from the MV model and from JIMWLK evolution. The magnitude and transverse 
momentum dependence of the v„{2} values suggest that the fluctuations present in the initial fields 
are a significant contribution to the observed anisotropies. 
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I. INTRODUCTION 

One of the most surprising results from the LHC 
proton-nucleus collision experiments have been the kind 
of azimuthal multiparticle correlation structures [1-7] 
(see also RHIC results from deuteron-gold collisions 
[8, 9]) that have, in larger collision systems, been at- 
tibuted to hydrodynamical flow. The particle multiplici¬ 
ties in these collisions systems are large enough for some 
collective effects to take place. Many of these structures 
have indeed been successfully descibed by hydrodynam¬ 
ical calculations [10, 11]. This agreement requires, how¬ 
ever, a very specific model of the geometry of the initial 
state [12]. It is also not clear whether these small sys¬ 
tems are within the regime of validity of a hydrodynami¬ 
cal description with realistic values of the energy density, 
viscosity and system size [13]. 

The primary collisions leading to energy deposition 
in the central rapidity region are, at the high energies 
reached at the LHC, characterized by very strong non¬ 
linear color fields [14]. These fields are, to leading order 
in the coupling constant, boost invariant. This immedi¬ 
ately leads to the presence of long range azimuthal corre¬ 
lations in particle production [15-22]. In larger collision 
systems, the structure of these correlations in azimuthal 
angle and transverse momentum is strongly influenced 
by collective behavior in the later evolution stages of the 
system. However, in smaller systems, such as proton- 
nucleus collisions, these collective effects are presumably 
less significant than in nucleus-nucleus collisions. This 
raises the intriguing possibility that in proton-nucleus 
collisions also the azimuthal structure of the initial stage 
color fluctuations could be directly visible in the measur¬ 
able particle spectrum. 

We do not yet have a very solid quantitative under¬ 
standing of the relative importance of initial color field 
and later evolution effects for generating anisotropies in 
particle production. A complete calculation of azimuthal 
anisotropies in this context requires complicated mod¬ 
eling that includes the color field and nucleonic scale 


fluctuations in the nucleus [23] and in the proton [24], 
combined with a calculation of the time evolution of the 
initial color fields and eventual matching to a hydrody¬ 
namical description [25]. We will not attempt to carry 
out this whole program here, but concentrate in this pa¬ 
per only on a part of it, namely the anisotropies produced 
when a bunch of valence quark-like particles in the fun¬ 
damental representation of the gauge group scatter off 
the color field of a large nucleus. The physical picture 
(see [20, 21] and more recently [26-29]) in our calcula¬ 
tion is that of valence quarks from the probe deflected in 
a preferred transverse direction by a domain in the target 
color field. This generates a multiparticle correlation that 
probes the spatial fluctuations of the target. Our calcula¬ 
tion extends the work in [28, 29] in two significant ways. 
Firstly, we perform the Fourier-transform from coordi¬ 
nate to momentum space, in order to get an azimuthal 
harmonic coefficient corresponding to real produced par¬ 
ticles. Secondly, unlike [28, 29], we correlate the particles 
in a given p^-bin with a reaction plane determined by 
all the produced particles using the 2-particle cumulant 
method. 


II. AZIMUTHAL CORRELATIONS IN CGC 
FIELDS 

In the “hybrid formalism” for particle production 
in the dilute-dense limit [30-33] the quark spectrum 
produced in a collision is proportional to the Fourier- 
transform of the two point function of Wilson lines in 
the color field of the target 

PC f p-Pr-(xT-yT) J_ Y ('1') 

XT,yT 

The Wilson lines V{xt) in (1) are, in the CGC descrip¬ 
tion, stochastic random SU(3) matrices in the represen¬ 
tation of the projectile. To calculate the single inclusive 



2 


cross section one has to average Eq. (1) by the appropri¬ 
ate probability distribution of Wilson lines. 

In the leading order CGC treatment that we use 
here, multiparticle correlations can be calculated from 
the higher order moments of the Wilson line operators 
in Eq. (1). This corresponds to the so called “glasma 
graphs” [18, 19, 22], which have a very clear interpreta¬ 
tion in the hybrid formalism. The target nucleus is rep¬ 
resented by a sheet of color magnetic and color electric 
fields, which have a domain structure with a character¬ 
istic length scale 1/Qs in the transverse plane. When a 
small enough probe (comparable in size to the domain 
size) hits this target, the resulting particle production is 
has a preferred direction given by the direction of the 
color field in the domain. Since this direction fluctuates 
from event to event, there is of course no anisotropy on 
average, but the existence of a preferred direction in in¬ 
dividual events shows up in a global angular correlation 
among all of the produced particles, similarly to hydro- 
dynamical flow. We are neglecting here “connected” or 
“BFKL”-like correlations [22], that give rise to a back- 
to-back peak in the two-particle correlation. These cor¬ 
relations are typical “nonflow” correlations that involve 
only a few particles, which the experimental analyses of 
azimuthal anisotropy try to exclude. We will not discuss 
them further here, see however Refs. [22, 34] for more 
studies on these lines. 

It is evident from the above discussion that we expect 
the correlation to be very sensitive to the transverse size 
of the probe. In the case of calculating the initial condi¬ 
tion for an ion-ion collision the probe is large, with the 
consequence that the correlation is washed away by the 
sum over many independent domains in the transverse 
plane. Thus, in contrast to the correlations generated 
by collective flow, the effect discussed here is stronger in 
small collision systems than in large ones. 

The purpose of this paper is to analyze the azimuthal 
correlation structure of particle production using Eq. (I) 
in more detail. In particular, we want to study its de¬ 
pendence on the harmonic n, transverse momentum, and 
the transverse size of the probe. The practical procedure 
used here is the following. We first divide the pt range 
accessible on the lattice into bins. We use here 50 bins, 
but we have checked that the results are independent of 
the size of the bin. We then define the Fourier coefficient 
of the single particle spectrum as 


|pT|Gbin xt.Yt 

X - hT)Sp{yT - br)^ TrEj^Ey,. (2) 

The transverse coordinate profile of the probe has been 
taken as a Gaussian 


I —(x-r — I 


( 3 ) 


around an impact parameter by chosen randomly in the 
transverse plane of the target. We will present results for 
different values of the parameter B characterizing the size 
of the probe. Note that the coefficients (2) need not be 
normalized, since we will eventually divide by the angular 
average spectrum bg to construct the Fourier harmonic 
coefficient. We want to calculate the angular correlations 
with respect to a an event plane defined by all the all the 
produced particles, which form the “reference” that we 
correlate individual particles with. This is done following 
the procedure used in the experimental analysis (see e.g. 
the 2-particle cumulant method in [4]). For this we need 
to calculate also the reference coefficients 


6„(ref) = J d^pre'^^PT 

XT,yT 

X Spi^T - bT)5p(yr - br) ^ Tr E^^Ey, (4) 

integrated over all momenta. 

The target Wilson lines are drawn from a completely 
homogenous and isotropic distribution that fills the whole 
transverse lattice with periodic boundary conditions, and 
the probe is azimuthally symmetric. Thus there is no ge¬ 
ometrical (i.e. originating in the shape of the probe or 
the target) origin for azimuthal anisotropy present in the 
calculation. Since the probability distribution of Wilson 
lines is azimuthally symmetric (although the individual 
configurations are not), the correlations among the coef¬ 
ficients bn are diagonal: 

{K{PT)bm{qT)) oc (5m,„, (5) 

where () denotes averaging over the configurations of Wil¬ 
son lines in the target. Note that the single particle spec¬ 
trum Eq. (1) is explicitly real, configuration by configu¬ 
ration, leading to This can be shown by taking 

the complex conjugate of Eq. (1) and exchanging the in¬ 
tegration variables xt and yr- The two particle pair 
correlation function is now 

£ {KiPT)bn{qT))cos{nAip). (6) 


J g-iPT-i^T-yr) 


From this we can identify the correlation function Fourier 
coefficients (using the notation of [4]) 


VnA{PT, qr) 


{bn{PT)bn{qT)) 
{boiPT)bo{qT)) ’ 


( 7 ) 


and define the 2-particle cumulant azimuthal harmonic 
as in [4] as 


Vn{2} 


VnA{PT,Ki) 

^E„A(ref,ref) 

(&n iPT)brviief)) 

{boiPT)bodef)) 


(&* (ref)6„(ref)) 
(b2(ref)6o(ref)) 


( 8 ) 


Sp{xT - hr) = exp 


( 9 ) 
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A nice interpretation of Eq. (9) can be obtained by writ¬ 
ing it as a product of three terms, 


Vn{2} = u„{bp} 


-^»(PT,ref) 

i?o(pT,ref) ■ 


Here we denote by 


^^n{bp}^ 


{K{pt)K{pt)) 

{hl{pT)ho{pT)) 


( 10 ) 


( 11 ) 


the flow coefficient for particles in the pt bin with re¬ 
spect to the event plane of that pt bin (“bp” stands for 
“bin plane”). This is the equivalent (although here in 
momentum, not position space) of the quantity calcu¬ 
lated in [28]. This is then corrected by two “correlation 
coefficients”. The first one is the correlation coefficient 
between the reference reaction plane and the pr-bin re¬ 
action plane: 


i?„(pT, ref) 


(^n(PT)&n(ref)) ^ ^ 

i.PT)K{ pt)) {bl (ref)&„(ref)) “ 


( 12 ) 

where the inequality follows from the Schwartz inequal¬ 
ity. The interpretation of this correction is clear: for 
a fixed anisotropy with respect to the pT-bin reaction 
plane, a decorrelation of the py-bin reaction plane from 
the reference reaction plane decreases the flow coefficient 
Vn{2}. The other correlation coefficient factor in (10) 


i?o(PT,ref) 


_ {b*o{PT)bo{ief)) _ ^ ^ 

\/(^S(?t)&o(pt)) (62(ref)6o(ref)) “ 


(13) 

is related to the multiplicity and appears in the denom¬ 
inator, increasing z;„{2}. This can be understood as fol¬ 
lows: with larger fluctuations in the pT’-bin multiplicity 
that are independent of the reference multiplicity, a hxed 
correlation between bnipr) and bn{Tei) implies a larger 
correlation between p^-bin and reference reaction planes. 
In other words, since ^ Vnbo, for a given correlation 
between (ref and px) b^s, the smaller the correlation be¬ 
tween bo’s, the larger must the correlation between u„:s 


be. 


We take the Wilson lines V{xt) appearing in Eq. (1) 
either from the MV model or resulting from JIMWLK 
evolution of the distribution of Wilson lines. Both are 
discretized on a 1024^ transverse lattice. For the MV 
model we use a (fundamental representation) saturation 
scale of Qstt = 0.119, where a is the lattice spacing. 
The JIMWLK calculation starts with a MV model at 
Qstt = 0.0220 and, after y = 10 units of evolution in ra¬ 
pidity (with running coupling) ends up with Qga = 0.117. 
The MV model Wilson lines are constructed following the 
procedure decribed in more detail in [35] and the running 
coupling JIMWLK evolution performed using the algo¬ 
rithm of [36] . The parameter values used here are exactly 
the same as for the 1024^-lattice in [37]. Note that we 
are only averaging two-point functions of the coefficients 
bn, not ratios ^ bn/bo- This makes the averaging proce¬ 
dure numerically quite stable and is physically the correct 


thing to do, since the pair correlation function Eq. (6) is 
the correct inclusive observable to be obtained via the 
target average [17, 38]. 


III. RESULTS AND DISCUSSION 

The numerical evaluations of the first four anisotropy 
coefficients are shown in Figs. 1, 2 and 3. The results are 
presented in scaling units, as plotted against pxlQs 
for different size probes, i.e. different VBQs- Also 
shown are the “bin plane” coefficients u„{bp}, defined 
by Eq. (11). To set the scale of the parameters in phys¬ 
ical units we note that the fundamental representation 
saturation scale around midrapidity at the LHC should 
be [39] around Qs ~ 1 GeV and the typical size of a 
proton in hard particle production at small x around 
H « 4 GeV~^ [40]. Thus a realistic probe size for LHC 
pA collisions would very roughly be VBQs « 2. 

The first immediate observation from the numerical re¬ 
sults is that the color field fluctuations indeed generate 
anisotropies that are large, of the order of the experimen¬ 
tally measured anisotropy coefficients. It seems therefore 
plausible that the color held huctuations do play a size¬ 
able role in the observed anisoptropy in small systems, 
and must be taken into account together with the how 
contribution. Also the momentum distribution has the 
same structure as the observed transverse momentum de¬ 
pendence of the how, hrst rising until ^ Qg and then 
decreasing. The “bin plane” coefficients u„{bp} do not 
decrease nearly as fast at high momentum, from which 
one can deduce that the decrease of the anisotropy co¬ 
efficients at large px follows from the decorrelation of 
the event plane in the px bin from the reference. This 
explains why this decrease was not seen in [28], where 
this decorrelation was not taken into account. The MV 
model has a gluon spectrum that is more sharply peaked 
around Qs, i.e. a narrower distribution of diherent size 
color held domains. This shows up in signihcantly larger 
values for the Vn coefficients. The main effect of JIMLK 
evolution is to add more small color held domains (larger 
Px gluons), which decrease the anisotropy of the particle 
spectrum. 

There is, however, an important caveat concerning any 
direct comparison of these results to experimental values. 
Namely, we were considering, in Eq. (1), only incoming 
quarks. For antiquarks one must replace the Wilson line 
by its Hermitian conjugate, which changes the sign of 
for odd n. Away from the very forward valence region 
in the probe, there are an approximately equal amount 
of quarks and antiquarks present, with contributions to 
that therefore cancel. Gluons do not have nonzero 
odd harmonics in this mechanism, because the adjoint 
representation is real and thus odd 6„’s vanish. Any odd 
harmonic surviving in the hnal state around midrapidity 
must therefore have an origin that is different from the 
one discussed here. 

The other word of caution in interpreting these results 
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FIG. 1: Second harmonic coefficient V2{2} calcnlated with JIMWLK-evolved (left) and MV-model (right) Wilson line config¬ 
urations. The thin lines represent the cofficients n 2 {bp} (see Eq. (11)) calculated with respect to the event plane in the pr bin 
only. 




FIG. 2: Third harmonic coefficient V3{2} calculated with JIMWLK-evolved (left) and MV-model (right) Wilson line configu¬ 
rations. The thin lines represent the cofficients U 3 {bp} (see Eq. (11)) calculated with respect to the event plane in the pr bin 
only. 




FIG. 3: Fourth harmonic coefficient V4{2} calculated with JIMWLK-evolved (left) and MV-model (right) Wilson line config¬ 
urations. The thin lines represent the cofficients U 4 {bp} (see Eq. (11)) calculated with respect to the event plane in the pr bin 
only. 
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is related to the dependence on the size of the probe, 
parametrized here by the width of the Gaussian ^/B. As 
anticipated, the magnitude of the correlation, and its de¬ 
pendence on the transverse momentum, depends strongly 
on the size of the interaction region. Although one can 
quite well estimate this, it depends on nonperturbative 
physics in the proton and cannot ultimately be controlled 
in a weak coupling calculation. 

Results for azimuthal correlations in a full GYM sim¬ 
ulation have also been presented recently by Schenke, 
Schlichting and Venugopalan [41, 42]. Their calculation 
includes effects of both color field and nucleonic fluctu¬ 
ations in the probe proton and the target nucleus. The 
probe and target geometries also have a significant ef¬ 
fect through the GYM pre-equilibrium version of the 
usual hydrodynamical mechanism that converts spatial 
anisotropy to momentum space, leading to also odd har¬ 
monics. These geometrical effects have not been in¬ 
cluded in our work, which should therefore not be com¬ 
pared directly with experimental data. Our focus here 
has been, in stead, on quantifying the generic obser¬ 
vation that fluctuating color fields result in azimuthal 
anisotropies in multiparticle correlations, even in the ab¬ 
sence of anisotropies in the impact parameter depen¬ 
dence. 

As a conclusion, we have here studied the momentum 
space azimuthal anisotropy structure of the “color glass” 


gluon fields in a high energy nucleus, as they are seen 
by a small probe consisting of valence-like quarks. We 
also quantified here the effect of correlating the parti¬ 
cles with the event plane determined by all the produced 
particles, using the two-particle cumulant method at the 
parton level. The quantitative results strongly depend on 
the details of the pT’-distribution of gluons in the CGC 
wavefunction and on the transverse size of the probe. 
However, all the results show large contributions to the 
harmonics from these purely initial state effects. For odd 
harmonics they largely cancel between quarks and anti¬ 
quarks, but for even harmonics these are sizeable effects 
that need to be considered when interpreting the exper¬ 
imental results from proton-nucleus collisions. 
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